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Abstract 

Cardiac alternans, a beat-to-beat alternation in action potential duration (at the cellular 
level) or in ECG morphology (at the whole heart level), is a marker of ventricular fibrillation, 
a fatal heart rhythm that kills hundreds of thousands of people in the US each year. Inves- 
tigating cardiac alternans may lead to a better understanding of the mechanisms of cardiac 
arrhythmias and eventually better algorithms for the prediction and prevention of such dreadful 
diseases. In paced cardiac tissue, alternans develops under increasingly shorter pacing period. 
Existing experimental and theoretical studies adopt the assumption that alternans in homoge- 
neous cardiac tissue is exclusively determined by the pacing period. In contrast, we find that, 
when calcium-driven alternans develops in cardiac fibers, it may take different spatiotemporal 
patterns depending on the pacing history. Because there coexist multiple alternans solutions for 
a given pacing period, the alternans pattern on a fiber becomes unpredictable. Using numerical 
simulation and theoretical analysis, we show that the coexistence of multiple alternans patterns 
is induced by the interaction between electrotonic coupling and an instability in calcium cycling. 

1 Introduction 

Sudden cardiac death, attributable to unexpected ventricular arrhythmias, is one of the leading 
causes of death in the US and kills over 300,000 Americans each year [TJ. The induction and main- 
tenance of ventricular arrhythmias has been linked to single-cell dynamics [2J [3]. In response to an 
electrical stimulus, cardiac cells fire an action potential [JJ, which consists of a rapid depolarization 
of the transmembrane voltage (V m ) followed by a much slower repolarization process before return- 
ing to the resting value (Fig. [1]). The time interval during which the voltage is elevated is called 
the action potential duration (APD). The time between the end of an APD and the beginning of 
the next one is called the diastolic interval (DI). The time interval between two consecutive stimuli 
is called the basic cycle length (BCL). When the pacing rate is slow, a periodic train of electrical 
stimuli produces a phase-locked steady-state response, where each stimulus gives rise to an identical 
action potential (1:1 pattern). When the pacing rate becomes sufficiently fast, the 1:1 pattern may 
be replaced by a 2:2 pattern, so-called electrical alternans (SJ [6], where the APD alternates between 
short and long values. Recent experiments have established a causal link between alternans and 
the risk for ventricular arrhythmias [3 El [9] . Therefore, understanding mechanism of alternans is 
a crucial step in detection and prevention of fatal arrhythmias. 
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Figure 1: Schematic action potential showing the response of the transmembrane voltage to periodic 
electrical stimuli. 

Cellular mechanisms of alternans have been much studied. Summaries on this topic can be 
found in recent review articles by Shiferaw et al. [TO] and Weiss et al. At the cellular level, 

cardiac dynamics involves bidirectional coupling between membrane voltage (V m ) dynamics and 
intracellular calcium (Ca- m ) cycling. During an action potential, the elevation of V m activates L- 
type Ca currents Ic& to invoke the elevation of [Caj n ], which in turn triggers Ca release from the 
sarcoplasmic reticulum (SR), a procedure known as calcium-induced-calcium release (CICR) [12] . 
The V m ^Caj n coupling satisfies graded release, where a larger DI leads to an increase in the Ca 
release at the following beat since it allows more time for L-type Ca channels to recover. On the 
other hand, Ca release from the SR affects the APD in two folds: to curtail the APD by enhancing 
the inactivation of L-type Ca currents Ica'i and to prolong the APD by intensifying Na + /Ca 2+ 
exchange currents /NaCa- Therefore, depending on the relative contributions of /ca and -^NaCai an 
increase in Ca release may either shorten the APD (negative Ca; n — >V m coupling) or lengthen 
the APD (positive Ca; n -^V m coupling) [131 [T4l I15j . There exist two main cellular mechanisms 
of alternans. Firstly, alternans may be attributed to steep APD restitution [5], which is due to 
a period-doubling instability in the V m dynamics. In this case, Ca; n transient alternans, as a 
slave variable, is induced because V m regulates [Cai n ] via the L-type Ca currents and the sodium- 
calcium exchange currents. Secondly, alternans may be caused by a period-doubling instability in 
Cai n cycling, which is associated with a steep relationship between the sarcoplasm reticulum (SR) 
release and SR load [TOJ [TTj. In this case, APD alternans is a secondary effect via Cai n ^V m 
coupling. For ease of reference, we call the first mechanism APD-driven alternans and the second 
Ca-driven alternans. Interestingly, APD and Ca; n transient alternans can be electromechanically 
(E/M) in phase or out of phase [18|. I19j . In E/M in-phase alternans, a long-short-long APD pattern 
corresponds to a large-small-large [Caj n ] pattern, See Fig. [2] (a). In contrast, in E/M out-of-phase 
alternans, a long-short-long APD pattern corresponds to a small-large-small [Ca; n ] pattern, See Fig. 
[2](b). When alternans happens in isolated cells, the bidirectional coupling between APD and Caj n 
transient determines the relative phase of APD and Cai n transient alternans. In particular, APD- 
driven alternans always leads to E/M in-phase alternans whereas Ca-driven alternans is E/M in 
phase for positive Caj n — >V m coupling and out of phase for negative Cai n — >V m coupling [T3 | [III [T5] . 

The mechanism of alternans in multicellular tissue is more complicated since it involves electro- 
tonic coupling and conduction velocity restitution. Of particular interest is a phenomenon called 
spatially discordant alternans, in which different regions of the tissue alternate out of phase. Discor- 
dant alternans is arrhythmogenic because it forms a dynamically heterogeneous substrate that may 
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Figure 2: Schematic illustration of electromechanically in-phase (a) and out-of-phase (b) alternans. 

promote wave break and reentry [2, 3j. To study the spatiotemporal patterns of alternans, Echebar- 
ria and Karma derived amplitude equations that are based on APD-driven alternans [201 [2Tj . These 
amplitude equations not only are capable of quantitative predictions but also provide insightful un- 
derstandings on the arrhythmogenic patterns. In a recent article, Dai and Schaeffer [22] analytically 
computed the linear spectrum of Echebarria and Karma's amplitude equations for the cases of small 
dispersion and long fibers. 

Spatial patterns of alternans have been investigated in experiments [7J EU [2U [25] . Recently, 
Aistrup et al. [26] used single-photon laser-scanning confocal microscopy to measure Ca signaling 
in individual myocytes. They found that Ca alternans is spatially synchronized at low pacing rates 
whereas dyssynchronous patterns, where a number of cells are out of phase with adjoining cells, arise 
when the pacing rate increases. Aistrup et al. also observed subcellular alternans at fast pacing, 
where Ca alternans is spatially dyssynchronous within a cell. Using simulations of 1-d homogeneous 
tissue, Sato et al. [15] found that, in cardiac fibers with negative Ca— >V m coupling, Ca alternans 
reverses phase over a length scale of one cell whereas, in fibers with positive Ca— >V m coupling, 
Ca alternans changes phase over a much larger length scale. They interpreted this difference by 
showing that negative Ca— >V m coupling tends to desynchronize two coupled cells while positive 
Ca— >V m coupling tends to synchronize the coupled cells. 

Motivated by the aforementioned experimental and theoretical work, this paper aims to ex- 
plore spatial patterns of cardiac alternans. Through extensive numerical simulations, we find that 
complex spatial patterns of Ca alternans with phase reversals in adjacent cells can happen in ho- 
mogeneous fibers with both negative and positive Ca— >V m couplings. Most surprisingly, we find 
that the spatiotemporal pattern of cardiac alternans is not determined by the pacing period alone. 
Specifically, when calcium-driven alternans develops in multicellular tissue, there coexist multiple 
spatiotemporal patterns of alternans regardless of the length of the fiber, the junctional diffusion 
of Ca, and the type of Ca— >V m coupling. We further investigate the mechanism that leads to the 
coexistence of multiple alternans solutions. Our analysis shows that multiple alternans solutions 
are induced because of the interaction between electrotonic coupling and an instability in Cai n 
cycling. 



2 Model Description 

2.1 Membrane dynamics 

We adopt a model of membrane dynamics that combines the calcium dynamics model developed 
by Shiferaw et al. [27] and the canine ionic model by Fox et al. [28]. In the following, we will 
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refer to this model as the Shiferaw-Fox model. Detailed formulations of the model can be found in 
[27\ 129] . The Shiferaw-Fox model has adopted two sets of parameters in the calcium dynamics to 
account for negative and positive Cai n — >V m couplings. Besides the phase difference, the two sets of 
parameters also produce alternans at different values of BCL. Using Shiferaw's default parameters, 
we find alternans happens at BCL~ 401 ms for negative Ca; n — >V m coupling and BCL~ 323 ms 
for positive Cai n — >V m coupling. Fig. [3] shows the bifurcation diagrams in APD and peak value of 
[Cai n ] for negative Ca; n ^V m coupling. The bifurcation diagrams for positive Cai n ^V m coupling 
are similar and thus are not shown here. We note that, in simulations of isolated cells using the 
Shiferaw-Fox model, alternans solutions do not depend on the initial condition nor on the pacing 
history. However, as we will show in the following, fibers based on the Shiferaw-Fox model possess 
multiple alternans solutions, which are sensitive to the initial condition and the pacing protocol. 
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Figure 3: Bifurcation diagrams of an isolated cell using the Shiferaw-Fox model. 



2.2 Simulation of fiber models 

We study paced, homogeneous fibers, which can be modeled using the cable equation: 

dv d 2 v 1 

m- D d^~cZ [Iion + Icxth (1) 

where v represents V m , D = 5 x 10 -4 cm 2 /ms represents the effective diffusion coefficient of v in the 
fiber, C m = 1 fiF/cm 2 represents the transmembrane capacitance, /i on is the total ionic current, 
and I ex t represents the external current stimulus. The ionic current Ji on is computed using the 
Shiferaw-Fox model. The current stimulus has duration 1 ms and amplitude 80 fxA/fxF. This paper 
studies fibers of various lengths. For two coupled cells, we pace the left cell. For longer fibers, we 
pace the leftmost few cells to ensure propagation. For example, the leftmost 5 cells are paced in 
simulating a fiber of 100 cells. The cable equation (TTJ) is solved using the finite difference method 
with a space step of Ax = 0.015 cm and time step of At = 0.1 ms. No-flux boundary conditions 
are imposed at both ends of the fiber 



2.3 Pacing protocols 

To study the onset and development of alternans, we pace both single cells and fibers of various 
lengths with several pacing protocols, which are briefly described below. 
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(i) In the downsweep protocol [30], the cell/fiber is paced periodically with period BCL until 
it reaches steady state. Then, the pacing period is reduced by AB and the procedure is repeated 
many times. Note that this protocol is also known as dynamic pacing protocol |31j . 

(ii) The perturbed downsweep protocol, proposed by Kalb et al. [30], can be regarded as a 
perturbation to the downsweep protocol. At each pacing period BCL, the cell/fiber is first paced 
N beats to reach steady state. Then, a longer pacing period is applied at the N+lst pacing, after 
which the original pacing period is applied for 10 beats to allow the tissue to recover its previous 
steady state. Next, a shorter pacing period is applied and followed by 10 beats of the original 
pacing period. Finally, the pacing period is reduced by AB and the procedure is repeated. 

(iii) To explore the possibility for multiple alternans solutions, we set up certain initial condition 
and pace the tissue with period BCL to reach steady state, a process we call direct pacing. 

(iv) To explore the origin of an alternans pattern, we use the upsweep protocol [30], which is a 
reversed downsweep protocol. 

3 Spatiotemporal Patterns of Alternans: Numerical Exploration 

We simulate fibers using the Shiferaw-Fox model with both negative and positive Caj n — ^V m cou- 
plings under various conditions. Default parameters in Shiferaw's code [29] are used unless other- 
wise specified. Despite quantitatively significant differences, we find both types of couplings lead 
to the coexistence of multiple alternans solutions. For clarity, we start with the results for negative 
Cai n ^V m coupling and defer the results for positive Caj n — >V m coupling in a later subsection. 

3.1 Coexistence of multiple solutions 

We first consider a homogeneous fiber of 100 cells with negative Caj n — >V m coupling. To our 
surprise, numerical simulations show that when the fiber is in alternans, there coexist multiple 
solutions for a given pacing period. For example, Fig. [J] shows 6 selected solutions of alternans 
for the fiber paced at BCL=375 ms. Here, the steady-state solutions in panels (a-c) are obtained 
using the downsweep protocol with step size AB=1 ms, 2 ms, and 25 ms, respectively. The pacing 
protocols are started from BCL=500 ms in (a) and (c) and from BCL=499 ms in (b). We note that 
the solution of a downsweep protocol is not influenced by the initial condition at the starting, long 
BCL; instead, the solution is sensitive to the step size AB. The steady-state solutions in panels 
(d-f) are obtained by pacing the fiber at BCL=375 ms with prescribed initial conditions for 200 
beats, the so-called direct pacing. The initial condition of the fiber in panel (d) is uniform, i.e., all 
cells are assigned the same resting voltage, gating variables, and ionic concentrations. The initial 
condition in panel (e) is same as that in (d) except that [Caj n ] is assigned to be 0.54 fjM. for the 
first 35 cells and 0.66 for the remaining cells. Interestingly, this initial condition leads to a 
steady-state [Caj n ] pattern, which, besides the phase reversal between cells 35 and 36, has another 
phase reversal between cells 11 and 12. The initial condition in panel (f) is the same as that in (d) 
except that [Cai n ] is randomly assigned for cells on the fiber according to a uniform distribution 
in the interval of 0.45 /iM to 0.75 //M. In all protocols, we pace the fiber for 200 beats at each 
BCL and plot the last 10 beats at BCL=375 ms. The simulation results in Fig. [5] demonstrate 
that the alternans on a fiber is not solely determined by the pacing period. Instead, the solution 
is sensitive to the pacing protocol and the initial condition. It is worth noting that, besides the 
solutions shown in Fig. |4j there exist many other solution patterns. In particular, there exist many 
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complex patterns similar to Fig. H] (f). In the following, we will verify whether the phenomenon is 
influenced by the length of the fiber, junctional Ca diffusion, or Cai n -^V m coupling. 
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Figure 4: Six selected alternans solutions for a homogeneous fiber of 100 cells with negative Ca^ V m 
coupling when paced at BCL=375 ms. Panels (a)-(c) are obtained using downsweep protocols from 
500 ms to 375 ms in steps of 1, 2, and 25 ms, respectively. Panels (d)-(f) are obtained by directly 
pacing at BCL=375 ms with different initial distributions in Cai n , see text for details. To verify if 
a solution reaches steady state, we plot the last 10 beats of the simulation results, where odd beats 
are represented by red solid lines and even beats by blue dashed lines. The coincidence of all odd 
beats and that of all even beats confirm that the solution is indeed in steady state. 



3.2 Influence of the length of the fiber 

We repeat the numerical simulations for fibers of various lengths and find that the phenomenon 
persists regardless of the length of the fiber. Probably, the most illuminating example is a "fiber" 
of two coupled cells. Denoting the voltage at cell 1 as V\ and that at cell 2 as V2, we can simulate 
the two coupled cells by integrating the following equations |32| 

(/ion + /cxt) + £ ( ^" 2 yi) , (2) 

Using various pacing protocols and initial conditions, we 
find the two coupled cells can have spatially desynchronized (Fig. [5] (a)) and synchronized (Fig. [5] 
(b)) alternans. Note that in the desynchronized pattern, APD in both cells exhibit a beat-to-beat 
variation of a few hundredth of ms, which is impossible to observe in experiments. 
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Figure 5: A "fiber" of two cells with negative Ca-> V m coupling can have both spatially desyn- 
chronized and synchronized alternans solutions. The pacing period is 375 ms. The last 10 beats of 
the steady-state solutions are plotted, where odd beats are represented by red solid lines and even 
beats by blue dashed lines. 



3.3 Influence of junctional Ca diffusion 

The original Shiferaw-Fox model does not include the diffusion of Ca between neighboring cells [33J . 
Physiologically, there exists gap junctional Ca diffusion although its magnitude is several orders 
lower than voltage diffusion [341 [35] and thus is typically neglected in cardiac modeling. One may 
wonder whether including junctional Ca diffusion will affect the simulation results. To answer this 
question, we introduce junctional Ca diffusion to the Shiferaw-Fox model by modifying the equation 
governing [Caj n ] as follows [32J: 

dt ~ c dx \dx rt Ui dx) c ' [ > 

where Cj represents [Cam], D c = 3 x 10~ 9 cm 2 /ms is the Ca diffusion coefficient [MIES], Z c is the 
valence of Ca, F is the Faraday constant, R is the gas constant, T = 300 K is the temperature, 
and I c represents the Ca; n currents in the Shiferaw-Fox model. With this modification, we repeat 
the numerical simulations and find there still exist multiple alternans solutions. 

For example, Fig. [6] shows three selected solutions for a homogeneous fiber of 100 cells paced at 
BCL=375 ms. Figures [6] (a) and (b) are obtained using downsweep protocols same as those in Fig. 
H](b) and (c), respectively. Figure [6] (c) is obtained using a random initial distribution in [Ca m ] (cf. 
Fig. 111(f)). Comparing results in Figs. 0] and [6] shows that the coexistence of multiple solutions is 
not influenced by junctional Ca diffusion. 



3.4 Influence of Ca in — >V m coupling 

Simulations show that fibers with positive Ca; n ^V m coupling also possess multiple alternans 
solutions. For example, Fig. [7] shows 3 selected solutions for a homogeneous fiber of 100 cells with 
positive Cai n — >V m coupling, paced at BCL=300 ms. Figure [7] (a) is obtained via a perturbed 
downsweep protocol, where the step size is AB=25 ms and a long and a short perturbations of 
±20 ms are applied at each BCL, see section 2 for details of the perturbed downsweep protocol. 
Figures [7] (b) and (c) are obtained by directly pacing the cell at BCL=300 ms using a uniform and 
a random initial distributions in [Cai n ], respectively. 
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Figure 6: Three selected alternans solutions for a homogeneous fiber model, where junctional Ca 
diffusion is included. The cells have negative Ca— > V m coupling and the pacing period is 375 ms. 
Panels (a) and (b) are obtained using downsweep protocols and panel (c) is obtained by direct 
pacing with random initial distribution of Caj n . The last 10 beats of the steady-state solutions are 
plotted, where odd beats are represented by red solid lines and even beats by blue dashed lines. 
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Figure 7: Three selected alternans solution in a homogeneous fiber of 100 cells with positive Ca— * V m 
coupling. The pacing period is 300 ms in all panels. Panel (a) is obtained using a perturbed 
downsweep protocol. Panels (b) and (c) are obtained by direct pacing with with initial uniform 
distribution and random distribution of Ca; n , respectively. (See text for details.) The last 10 beats 
of the steady-state solutions are plotted: odd beats are red solid and even beats are blue dashed. 



We simulate two coupled cells with positive Cai n — > V m coupling and find that alternans can 
be both spatially synchronized and spatially desynchronized, see Fig. [8l In Figs. [5]and[U values of 
APD are almost identical in the two cells because conduction time across a cell's length is negligible. 
Thus, these examples indicate that the coexistence of multiple solutions does not depend on steep 
CV restitution. We note that, using a downsweep protocol with step size of 2 ms, Sato et al. |15] 
also studied alternans in two coupled cells. They observed the desynchronized solution for negative 
Cai n -^Vj] coupling and the synchronized solution for positive Ca; n — >V m coupling. Using numerical 
simulation, Sato et al. showed that negative Cai n — >V m coupling tends to desynchronize two coupled 
cells whereas positive Cai n — >V m coupling tends to synchronize them. Thus, we hypothesize that 
the synchronized solution for negative Caj n — >V m coupling and the desynchronized one for positive 
coupling are induced by electrotonic coupling, which will be verified in the next section. 




Figure 8: (Color online) A "fiber" of two cells with positive Ca— > V m coupling can have both 
spatially synchronized and desynchronized alternans solutions. The pacing period is BCL= 300 ms. 
The last 10 beats of the steady-state solutions are plotted: odd beats are red solid and even beats 
are blue dashed. 



4 Origin of the Indeterminacy 

4.1 Numerically tracing an alternans solution 

The coexistence of multiple alternans solutions in a fiber is surprising because the underlying cell's 
model possesses a supercritical period-doubling bifurcation, describing a transition from a unique 
1:1 solution to a unique alternans solution (Fig. [3]). The transition to alternans in a fiber appears 
to be more complicated. First, when paced at sufficiently large values of BCL, a fiber has a single 
1:1 solution, regardless of the initial condition. However, as shown in the previous section, the 
fiber can develop multiple patterns of alternans. Therefore, it is interesting to ask how different 
alternans patterns are related to one another and how they are connected to the 1:1 solution. To 
address these questions, we utilize the upsweep protocol to trace the "origins" of different alternans 
solutions. Starting at BCL=375 ms from a pattern in Fig. HI we increase BCL by 1 ms every 
200 beats until the alternans solution either changes to a different alternans pattern or becomes a 
1:1 pattern. Then, we mark that value of BCL as the "starting" point for the studied alternans 
pattern. Under this protocol, each solution in Fig. H] transforms into another alternans solution at 
a starting BCL shown in Table [TJ Specifically, pattern (a) "starts" at BCL=416 ms and patterns 
(b)-(f) start at BCL values between 402 and 405 ms. Further numerical simulations show that the 



9 



fiber first undergoes alternans at BCL=425 ms. This is counterintuitive because, in the underlying 
cell's model, the onset of alternans occurs at 401 ms (cf. Fig. [3]). Since alternans solutions in the 
fiber can occur earlier than alternans in the single-cell model, we hypothesize that these solutions 
are induced from the interaction between cellular dynamics and electrotonic coupling (induced 
alternans). We further hypothesize that alternans may also result from the intrinsic instability 
mechanism that leads to alternans in the underlying cellular model (intrinsic alternans). Noticing 
that the spatially concordant alternans in pattern (d) starts at BCL=402 ms, near the onset of 
alternans in the single-cell model, we hypothesize the concordant alternans is a result of intrinsic 
alternans. 



Table 1 : Origins of the alternans patterns in Fig. [H 



pattern 


(a) 


(b) 


(c) 


(d) 


(e) 


(f) 


starting BCL (ms) 


416 


405 


404 


402 


403 


404 



Using upsweep protocols, we also trace the origins of the alternans in Fig. [5] for the two coupled 
cells. Numerical simulations show that the desynchronized solution in Fig. 0(a) stops at BCL=405 
ms and the synchronized solution in Fig. [5] (b) stops at BCL=401 ms. Moreover, we also find that 
different alternans solutions for fibers with positive Caj n ^V m coupling stop at different values of 
BCL. Thus, the bifurcation structure for a fiber is much more complicated than that for a single cell. 
Moreover, the numerical results suggest that the multiple alternans patterns born from different 
bifurcations are due to the interaction between electrotonic coupling and instability in calcium 
cycling. 

4.2 Bifurcation Mechanism 

Shiferaw et al. showed that the bidirectional V m /Ca coupling in the Shiferaw-Fox can be captured 
by a 2-d mapping model |13j . The work of Shiferaw et al. provides a theoretical framework 
for understanding various experimental observations including electromechanically in-phase/out- 
of-phase alternans and quasi-periodic oscillations of voltage and calcium. Here, we follow their 
approach and take a mapping model in the following form [36], 

A n+l = f (D n ) + 7 (C n+1 - C crit ) , (5) 
C n+ i = fiD n - g (C n ) , (6) 

where A n , D n , and C n represent the APD, DI, and peak [Ca; n ] concentration of nth beat, respec- 
tively. By definition, it follows that D n =BCL— A n , see Fig. [TJ Function / stands for the APD 
resitution. The second term in Eq. ([5]) accounts for the influence of [Ca; n ] on APD. Negative 
Cai n -^V m coupling corresponds to 7 < and positive coupling to 7 > 0. Due to graded release, 
H has to be positive [15J. Function g determines the relation between Ca concentrations in two 
consecutive beats. 

Using map ([5l [6]) , we will show how the interaction between electrotonic coupling and instability 
in calcium cycling leads to multiple alternans patterns in fibers. Since the analysis below does not 
depend on the forms of / and g, we do not specify their forms here. Instead we present numerical 
examples based on concrete forms of / and g in the Appendix. The numerical examples there 
show that fiber models based on the map (0 [6]) are able to reproduce the phenomena observed in 
simulations of the Shiferaw-Fox model. 
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4.2.1 Bifurcation for single cells 

Using a mapping model similar to (|5I6|) . Shiferaw et al. [T3] carried out a bifurcation analysis to 
show how different alternans solutions in single cells arise as a result of V m /Ca coupling. Here, we 
briefly review that analysis, which serves as a starting point to understand bifurcation for fibers. 
To this end, consider a paced cell described by map (|5l6p . A 1:1 solution is a fixed point of the 
map and its stability is determined by the Jacobian matrix: 

where /' measures the slope of the APD restitution and g' measures the slope of the [Ca; n ] relation, 
and all derivatives are evaluated at the fixed point. 

The strength of V m /Ca coupling varies among species |13| . In this paper, we assume a weak 
coupling between V m dynamics and Ca cycling, which allows easy analysis and yields insight on 
the instability mechanisms in cardiac cells as well as in fibers. To this end, we assume 7 fx « /', 
g' in the following. Using a perturbation technique [37], we find the eigenvalues of J to first order 
in 7 fi are 

Since 7 /i is small, one can see that a period-doubling bifurcation occurs if one of the slopes, /' 
or g' , becomes sufficient large (compared to 1). Thus, if the bifurcation is caused by increasing in 
/', we say the alternans is APD driven. On the other hand, it is called a Ca-driven alternans if 
the bifurcation is due to increasing in g' . For simplicity, if alternans is APD driven, we assume g' 
remains small and will not induce dynamic instability for all physiologically interesting parameter 
regimes; and vice versa. 

APD-driven alternans: In case of APD-driven alternans, it follows that /' > g' . To first 
order in 7 fi, the unstable eigenvector is 



1 _ T£f 
f'~9' {f'-g'f 



(9) 



which is in phase in APD and [Cai n ]. Therefore, APD-driven alternans is electromechanically in 
phase. 

Calcium-driven alternans: In case of Ca-driven alternans, it follows that g' > /'. To first 
order in 7 fi, the unstable eigenvector is 

I)"- (10) 

For positive Ca; n — >V m coupling, it follows that 7 > and thus alternans is electromechanically 
in phase. On the other hand, for negative Ca; n ^V m coupling, it follows that 7 < and thus 
alternans is electromechanically out of phase. 

4.2.2 Bifurcation for two coupled cells 

To understand the bifurcation mechanism for fibers, we start with the extreme case of two coupled 
cells. Denoting the locations of the cells by 37 and X2, we represent the APD, DI, and [Cai n ] of 
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cell % for the nth beat as A n (xj), D n (xj), and C n (xj), respectively. Following Fox et al. [38], we 
account for electrotonic coupling between the cells using weighted averaging: 

A n+1 On) = p (0) F (xi) +p(l)F (x 2 ) (11) 

A n+1 (x 2 ) = p (1) F (xi) + p (0) F (x 2 ) (12) 

C n+ i(xi) = pD n (x 1 ) -g(C n (x 1 )) (13) 

C n+ i(x 2 ) = pD n (x 2 ) -g(C n (x 2 )), (14) 

where 

F (xi) = / (D n fa)) + 7 (Cn+i (a*) - C crit ) , (15) 

and p (0) and p (1) are the weighting functions. Here, we have neglected the effect of dispersion due 
to conduction velocity restitution since the "fiber" is short and thus the influence of conduction 
restitution to APD is insignificant. On the time scale of one APD, V m diffuses on a spatial scale 
of a few tens of cells |20j. Thus, the coupling in APD between two neighboring cells is strong. As 
a result, p(0) « p(l) and p(0) +p(l) = 1. Therefore, we let p(Q) = 1/2 + e and p(l) = 1/2 - e, 
where < e « 1. For a 1:1 solution, the two cells have the same values of APD and [Cai n ]. Thus, 
the Jacobian matrix can be written as 



J 



/ p(0) (-/'-7M) p(l)(-/'-7«) -P(0)jg' "P(l)75'\ 

p(l)(-f'- 7 p) p (0) (-f-7/x) -p(l)7ff' -p(o)jg' 
-p 
V -u - 5 ' / 



(16) 



where all derivatives are evaluated at the 1:1 solution. 

Using perturbation theory [37], we find the eigenvalues to first order in 7 p and e are 



~f ( 1 + 7^7)' (i + T 2 ^)' -2e/', and - 5 ' -2e 7 ^ (17) 

Note that the first two eigenvalues are the same as those for the case of single cells (cf. Eq. [8]). 
Thus, they are due to intrinsic membrane dynamics. The last two eigenvalues are induced by the 
electrotonic coupling as manifested by the small coupling parameter e. 

APD-driven alternans: For APD-driven alternans, the intrinsic bifurcation is the only mech- 
anism to produce alternans. The corresponding eigenvector to first order in 7 u and e is 

Thus, APD-driven alternans tends to be electro-mechanically in phase and spatially synchronized. 

Calcium-driven alternans: For Ca-driven alternans, bifurcation can happen intrinsically 
when g'+ 7 p g' j {g' — /') = 1 {intrinsic bifurcation). Bifurcation can also be induced by electrotonic 
coupling when g' + 2^7 p = 1 (induced bifurcation). To first order in jp and e, the eigenvector 
corresponding to the intrinsic bifurcation is 

<?'-/" 9'-/" X ' 1 J ■ {Ly) 
To first order in 7 p and e, the eigenvector corresponding to the induced bifurcation is 

( 2e 7 a, -2e 7 p, 1, -if. (20) 
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Therefore, while the intrinsic bifurcation gives birth to a spatially synchronized solution, the in- 
duced bifurcation leads to a spatially desynchronized pattern. And the electromechanical phase is 
determined by the sign of 7. 

Subtracting the intrinsic eigenvalue from the induced one yields 7/i (g' (g' — /') — 2e), which 
has the same sign as 7 since g' > f (Ca-driven alternans) and [i > (graded release). If 7 < 0, 
the induced eigenvalue is more negative; therefore, the spatially desynchronized pattern, born from 
the induced bifurcation, occurs before the spatially synchronized pattern, born from the intrinsic 
bifurcation. On the other hand, if 7 > 0, the intrinsic eigenvalue is more negative; therefore, 
the spatially desynchronized pattern occurs after the spatially synchronized pattern. The analysis 
is in agreement with numerical simulations of the Shiferaw-Fox model using upsweep protocols. 
Therefore, when BCL is continuously decreased from a 1:1 solution, it is more likely to observe the 
spatially synchronized solution for the positive Caj n — >V m coupling case whereas it is more likely to 
observe spatially desynchronized pattern in the negative Cai n — >V m coupling case. This is probably 
why Sato et al. observed different solution patterns for negative and positive Ca; n — >V m couplings 
in their downsweep simulations. 

The conclusion based on two coupled cells may be extended to short fibers; however, we expect 
long fibers have more complicated phenomena than presented here. In particular, because the 
dispersion effect may play an important role in long fibers as shown by Echebarria and Karma 
[20\ |2T] and by Dai and Schaeffer [22] . Nevertheless, the simple case sheds light on the effect 
of electrotonic coupling, the difference between APD and Ca-driven alternans, and the difference 
between negative and positive Caj n ^V m couplings. 

4.2.3 Understanding bifurcation for fibers 

For Ca-driven alternans in fibers, we start with the extreme case of fi = 0. As can be seen from the 
cellular model (|5|6p . in this case, Caj n dynamics becomes independent of V m dynamics. Because 
Cai n dynamics is not influenced by electrotonic coupling, when alternans develops in a fiber, cells 
can arbitrarily choose their phases in [Caj n ]. Depending on the initial conditions, Ca alternans on 
a fiber can take different spatial patterns. When fi is small but nonzero, [Caj n ] in different cells are 
weakly coupled due to feedback of V m dynamics and electrotonic coupling. If this weak coupling 
does not suppress the coexistence of multiple solutions, the alternans pattern will become sensitive 
to initial conditions and pacing protocols. Under certain pacing protocols, Cai n alternans on a fiber 
may possess a complex pattern with multiple phase reversals, as shown in Figs. H] (f), [6] (c), and 
[7] (c). In these examples, the spatial patterns of APD are much less complex, which is because 
APD on the fiber is an averaging effect over many cells due to the fast diffusion of V m [20, 21J. 
As a result of electrotonic coupling, different spatial patterns of Caj n alternans may correspond 
to similar spatial patterns of APD. Moreover, as shown in simulations of the Shiferaw-Fox model 
(see Figs. [U El and ED and of the mapping model (see Fig. [HI] in Appendix), spatial dyssynchrony 
in Cai n alternans may induce spatially discordant APD alternans, which verifies the hypothesis of 
Aistrup et al. [26] . 

5 Summary and Discussion 

Using numerical simulation and theoretical analysis, we have investigated spatiotemporal patterns 
of calcium-driven alternans. The main finding is that although alternans of an isolated cell is 
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solely determined by the pacing period, the solution in a fiber is sensitive to pacing protocols and 
initial conditions. To the author's knowledge, this is the first report on the coexistence of multiple 
alternans solutions for cardiac fibers. We have further verified that the coexistence of multiple 
alternans solutions is independent of the length of the fiber, junctional Ca diffusion, and the type 
of Cai n - coupling. Since multiple solutions also exist for fibers of as few as a couple of cells, the 
phenomenon does not require steep conduction velocity restitution. Another interesting observation 
is that complex patterns of Caj n alternans with multiple phase reversals in neighboring cells may 
arise in a homogeneous fiber with both negative and positive Ca; n ^V m couplings. The simulation 
results also verify the hypothesis of Aistrup et al. that spatially desynchronized Ca signaling may 
lead to spatially discordant APD alternans [26]. 

We have also explored the bifurcation mechanism for the coexistence of multiple solutions. 
Numerical simulations and symbolic bifurcation analyses trace the onset of multiple alternans 
patterns to a number of bifurcations induced by the interaction between electrotonic coupling and 
an instability in calcium cycling. The bifurcation here bears similiarity to that of a ring model. 
For example, Courtmanche et al. [39] used the Beeler-Reuter model [40J to show that, when the 
length of the ring is reduced, there is an infinite-dimensional Hopf bifurcation, which gives birth 
to an infinite number of quasi-periodic solutions. Since the Beeler-Reuter model does not account 
for intracellular calcium cycling, it will be interesting to investigate how calcium-induced alternans 
propagates in a ring. 

Simulations in this paper lead to another interesting observation — the onset of altrenans in 
a fiber does not occur at the same value of BCL as that in the underlying single-cell model, a 
phenomenon against the common wisdom. This counter-intuitive observation indicates that one 
can not directly relate behavior of fibers to that of single cell models. Another intriguing example 
about the differences between alternans in fibers and alternans in single cells is presented in the 
work of Cherry and Fenton [H], who analyzed two models of canine ventricular myocytes: the Fox- 
MchHarg-Gilmour (FMG) model [28] and the Hund-Rudy (HR) model [12]. Simulations of both 
models show that the bifurcation structures for fibers are different than those for the corresponding 
single cells (see Fig. 5 in [41J). Most interestingly, the HR Od model (a single cell) shows alternans 
for BCL between 180 and 230 ms whereas the HR Id model (a 1.25 cm-long fiber) does not undergo 
alternans for all values of BCL studied. Examples in [H] as well as examples in the current paper 
show that stability of a fiber may not be inferred from stability of the single cell model and vice 
versa. 

Results in this paper are based on numerical simulations and hypothesized mapping models. 
Here, we describe potential experiments to verify the numerical observations. Probably, the simplest 
experiment would measure the spatial patterns of APD or Ca alternans on a fiber under various 
pacing protocols. It is well know that in vitro experiments suffer nonstationary drifts due to tissue 
dying. However, this is a slower process compared to typical experimental protocols, which last 
a few tens of minutes. Moreover, time constant of cardiac tissue is on the order of a few tens 
of seconds |30j . Therefore, it is possible to conduct multiple well-designed protocols in a time 
interval when tissue dynamics remains stationary. Other related experiments include the work of 
Fenton [13], where action potential is measured using microelectrode at a cell in a slice of paced 
dog epicardium tissue. Fenton showed that the onset and the form of alternans at the measured 
cell depend on the pacing history. While Fenton's observations may be attributable to factors such 
as cardiac memory, the coexistence of multiple alternans patterns in cardiac tissue may be another 
contributing element. 
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As pointed out by Qu and Weiss [H], because Caj n instability and V m instability are always 
being regulated by each other, it remains a challenge to identify the sources of instabilities in 
cardiac experiments. Recent work by Sato et al. [15] and by Jordan and Christini have 
developed theoretical criteria to assess the relative contributions of V m and Cai n dynamics in 
inducing cardiac alternans. The coexistence of multiple alternans patterns observed here may 
provide another possible criterion for this purpose. If the phenomenon was verified in experiments, 
it will raise questions in detection of alternans. New theories will also be needed to understand the 
development of different alternans patterns to ensure better prediction and control of spatiotemporal 
alternans. 
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A Simulation of the Mapping Model 

Here, we present numerical examples for single cells and for fibers using the mapping model (El EJ) - 
For simulation purpose, we let / (D n ) = a (l — /3 e~ D "/ T ) and g (C n ) = —s (C n — C cr it) — C m i n . 
Inspired by the results of Shiferaw et al. (Fig. 10 (d) in |27J), we choose the slope s in g to 
be s = | (1 — 0) s\ + ~ (1 + 6) S2, where 9 = tanh (k (C — C cr it) )• To produce numerical results, 
we adopt the following set of parameters: a = 500 ms, (3 = 0.62, r = 500 ms, 7 = —30 ms//xM, 
H = 0.0025 /xM/ms, C cri t = 1.6 fM, C min = 1.25 fM, si = -2.5, s 2 = 0.01, and k = 10 fiM" 1 . 
Note that we have arbitrarily chosen a negative value of 7 although a positive value of 7 will 
produce similar phenomena. Moreover, we make the magnitude of 7 fx small to simulate weak 
V m /Cai n coupling. 

A.l Alternans for single cells 

With the chosen parameters, map ([5j [6]) gives birth to alternans through a period-doubling bifur- 
cation at B = 410.958 ms, see the bifurcation diagrams in Fig. [9l For this simple model, it can 
be analytically shown that, for a given BCL, the cell has either a unique 1:1 solution or a unique 
alternans solution. However, as we will see below, a fiber has multiple alternans solutions induced 
from electrotonic coupling. 

A. 2 Alternans for fibers 

Based on the mapping model (J5j[6|) for single cells, we construct a coupled-maps model to simulate 
fibers. Coupled-maps models have been used by a few authors to study propagation of action 
potential [38 1 147 } UH l 149] . Here, we follow that approach. Specifically, we consider a fiber consisting 
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of M identical cells, located at Xi, % = 1, 2, . . . M. The distance between two neighboring cells is 
denoted by Ax. At cell Xi, APD and DI are related by the following equation: 

D n+ i (x^ = T n+ i (x^ - A n+ i (x^ , (21) 

where T n+ \ {xi) is the time interval between two consecutive activations of site X{. The time 
T n+ i (x^ is determined by the propagation time from the pacing site to Xi, that is, 

w^ + E^-E^, (22) 

where V n (xi) stands for the conduction velocity. We adopt a conduction velocity from Fox et 
al. [48]: V n (xt) = V max (1 - exp (- (D n (x^) + (3) /5)) with parameters V max = 0.72 cm/ms, /3 = 
17.408, and 5 = 14. To account for the effect of electrotonic coupling, we modify the weighted 
averaging formula in Fox et al. [38] as follows: 

^min(30, M-i) 



E7 =max'(-30,l- t ) W gj ^+1 (Xj) 
v min(30, M-i) , 
Z^j=max(-30,l-i) ; 



^n+U^j - ^ min (30, M-i) ' ^ 6 > 



where A n+ i (xj) = f (D n (xj)) + 7 (C n+ i (xj) — C cr jf) and w (j) = exp (— 0.0067 j 2 ). The equation 
for [Cam] reads 

C„+i = n D n ( Xi ) + 5 (C n ( Xi )) . (24) 

Simulations of the coupled-maps model (|21H24p show that the alternans pattern on a fiber 
depends on the pacing history, a phenomenon consistent with that observed in simulations of the 
Shiferaw-Fox model. For example, Fig. [10] shows 3 selected solutions for a homogeneous fiber of 
100 cells paced at B=375 ms. Panels in Fig. [10] have different initial conditions in [Ca; n ]. Panel (a) 
starts from a uniform initial distribution. Panel (b) starts from an initial condition, where [Ca; n ] 
is set to be 1.5 in the first 35 cells and 1.7 /xM in the remaining cells. Panel (c) starts from a 
random initial distribution. Thus, the coupled maps model is able to reproduce the coexistence of 
multiple alternans solutions as observed in the Shiferaw-Fox model. 
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Figure 10: Simulation of the coupled-maps model: three selected alternans solutions for a homo- 
geneous fiber of 100 cells with negative Ca— > V m coupling. The pacing period is 375 ms in all 
panels. Panel (a) starts from a uniform initial distribution in Ca; n . Panel (b) starts from an initial 
condition, where Cai n is set to be 1.5 fiM in the first 35 cells and 1.7 //M in the remaining cells. 
Panel (c) starts from a random initial distribution in Ca. The last 10 beats of the steady-state 
solutions are plotted, where odd beats are represented by red solid lines and even beats by blue 
dashed lines. 



References 



[1] |http : / / www . americanheart . org/, 

[2] A. Karma, 'Electrical alternans and spiral wave breakup in cardiac tissue', Chaos, 4, 461-472, 
1994. 

[3] A. Garfinkel, Y.H. Kim, O. Voroshilovsky, Z. Qu, J.R. Kil, M.H. Lee, H.S. Karagueuzian, J.N. 
Weiss, and P.S. Chen, 'Preventing ventricular fibrillation by flattening cardiac restitution,' 
Proc. Natl Acad Sci USA, 97:6061-6066, 2000. 

[4] R. Plonsey and R.C. Barr, Bioelectricity: A Quantitative Approach, Kluwer, New York, 2000. 

[5] J.B. Nolasco and R.W. Dahlen, 'A Graphic Method for the Study of Alternation in Cardiac 
Action Potentials,' Journal of Applied Physiology, 25:191-196, (1968). 

[6] A. Panfilov, 'Spiral breakup as a model of ventricular fibrillation,' Chaos 8:57-64, 1998. 

[7] J.M Pastore, S.D. Girouard, K.R. Laurita, F.G. Akar, D.S. Rosenbaum, 'Mechanism linking 
T-wave alternans to the genesis of cardiac fibrillation,' Circulation. 99:1385-1394, 1999. 

[8] D. M. Bloomfield, J. T. Bigger, R. C. Steinman, P. B. Namerow, M. K. Parides, A. B. Curtis, 
E. S. Kaufman, J. M. Davidenko, T. S. Shinn and J. M. Fontaine, 'Microvolt T-wave alternans 
and the risk of death or sustained ventricular arrhythmias in patients with left ventricular 
dysfunction,' J. Am. Coll. Cardiol., 47:456-463, 2006. 

[9] V. Shusterman, A. Goldberg and B. London, 'Upsurge in T-wave alternans and nonalternating 
repolarization instability precedes spontaneous initiation of ventricular tachyarrhythmias in 
human,' Circulation, 113:2880-2887, 2006. 



17 



[10] Y. Shifcraw, Z. qu, A. Garfinkel, A. Karma, and J.N. Weiss, 'Nonlinear Dynamics of Paced 
Cardiac Cells,' Ann. N.Y. Acad. Sci. 1080:376-394, 2006. 

[11] J.N. Weiss, A. Karma, Y. Shiferaw, P-S, Chen, A. Garfinkel, and Z. Qu, 'From Pulsus to 
Pulseless: The Saga of Cardiac Alternans,' 98:1244-1253, 2006. 

[12] M. Endo, 'Calcium release from the sarcoplasmic reticulum,' Physiol Rev 57:71-108, 1977. 

[13] Y. Shifcraw, D. Sato, and A. Karma, 'Coupled dynamics of voltage and calcium in paced 
cardiac cells,' Phys Rev E, 71:021903, 2005. 

[14] Y. Shiferaw and A. Karma, 'Turing Instability mediated by voltage and calcium diffusion in 
paced cardiac cells,' PNAS, 103:5670-5675, 2006. 

[15] D. Sato, Y. Shiferaw, A. Garfinkel, J.N. Weiss, Z. Qu, and A. Karma, 'Spatially Discordant 
Alternans in Cardiac Tissue: Role of Calcium Cycling,' Circulation Research, 99:520-527, 
2006. 

[16] E. Chudin, J. Goldhaber, A. Garfinkel, J. Weiss, and B. Kogan, 'Intracellular Ca(2+) dynamics 
and the stability of ventricular tachycardia,' Biophysics Journal, 77:2930-2941, 1999. 

[17] M.E. Diaz, S.C. O'Neill and D.A. Eisner, 'Sarcoplasmic Reticulum Calcium Content Fluctua- 
tion Is the Key to Cardiac,' Circulation Research, 94:650-656, 2004. 

[18] D.S. Rubenstein and S.L. Lipsius, Circulation 91, 201 (1995). 

[19] M.L. Walker and D.S. Rosenbaum, Cardiovasc. Res. 57, 599 (2003). 

[20] B. Echebarria and A. Karma, 'Instability and Spatiotemporal Dynamics of Alternans in Paced 
Cardiac Tissue,' Phys. Rev. Lett. 88:208101, 2002. 

[21] B. Echebarria and A. Karma, 'Amplitude equation approach to spatiotemporal dynamics of 
cardiac alternans,' to appear in PRE. 

[22] S. Dai and D. G. Schaeffer, 'Spectrum of a linearized amplitude equation for alternans in a 
cardiac fiber,' in process. 

[23] B.R. Choi, G. Salama, 'Simultaneous maps of optical action potentials and calcium transients 
in guinea-pig hearts: mechanisms underlying concordant alternans,' J Physiol. 529:171-188, 
2000. 

[24] E.J. Pruvot, R.P. Katra, D.S Rosenbaum, K.R. Laurita, 'Role of calcium cycling versus resti- 
tution in the mechanism of repolarization alternan,' Circ Res. 94:1083-1090, 2004. 

[25] R.P. Katra, E. Pruvot, K.R. Laurita, 'Intracellular calcium handling heterogeneities in intact 
guinea pig hearts,' Am J Physiol Heart Circ Physiol, 286:H648-H656, 2004. 

[26] G.L. Aistrup, J.E. Kelly, S. Kapur, M. Kowalczyk, I. Sysman-Wolpin, A.H. Kadish, and J. A. 
Wasserstrom, 'Pacing-induced Heterogeneities in Intracellular Ca Signaling, Cardiac Alter- 
nans, and Ventricular Arrhythmias in Intact Rat Heart,' Circulation Research, 99:65-73, 2006. 



18 



[27] Y. Shiferaw, M. Watanabe, A. Garfinkel, J. Weiss, A. Karma, 'Model of intracellular calcium 
cycling in ventricular myocytes,' Biophys J., 85:3666-3686, 2003. 

[28] J. J. Fox, J. L. McHarg, and R. F. Gilmour, Am. J. Physiol., 'Ionic mechanism of electrical 
alternans,' 282:H516-H530, 2002. 

[29] A Fortran code of the Shiferaw-Fox model was made by Yohannes Shiferaw and is available at 
http : //www. csun.edu/~yshiferaw/ shiferaw.html. 

[30] S.S. Kalb, H.M. Dobrovolny, E.G. Tolkacheva, S.F. Idriss, W. Krassowska and D.J. Gauthier, 
'The restitution portrait: a new method for investigating rate-dependent restitution,' J. Car- 
diovasc. Electrophysiol., 15:698-7. 

[31] M.L. Riccio, M.L., Koller, R.F. Gilmour Jr., 'Electrical restitution and spatio-temporal orga- 
nization during ventricular fibrillation,' Circulation Research, 84:955-963, 1999. 

[32] W. Krassowska, unpublished. 

[33] In models without junctional diffusion, [Ca; n ] in neighboring cells are indirectly coupled 
through their dependence on the voltage. 

[34] D.M. Bers, Excitation- Contraction Coupling and Cardiac Contractile Force (Developments in 
Cardiovascular Medicine), Kluwer, Boston, 2001. 

[35] N.L. Allbritton; T. Meyer, and L. Stryer, 'Range of Messenger Action of Calcium Ion and 
Inositol 1,4,5-Trisphosphate,' Science, 258:1812-1815, 1992. 

[36] Yohannes Shiferaw, Daisuke Sato, and Alain Karma, private communication. 

[37] A.H. Nayfeh, Perturbation Methods, John Wiley Interscience, New York, 1973. 

[38] J.J. Fox, M.L. Riccio, F. Hua, E. Bodenschatz, and R.F. Gilmour Jr., 'Spatiotemporal Tran- 
sition to Conduction Block in Canine Ventricle,' Circulation Research, 90:289-296, 2002. 

[39] M. Courtemanche, L. Glass, and J. P. Keener, 'Instabilities of a Propagating Pulse in a Ring 
of Excitable Media,' Physical Review Letters, 70:2182-2185, 1993. 

[40] GW. Beeler and H. Reuter, Reconstruction of the Action Potential of Ventricular Myocardial 
Fibres, J. Physiol. 268:177-210, 1977. 

[41] E.M. Cherry and F.H. Fenton, 'A Tale of Two Dogs: Analyzing Two Models of Canine Ven- 
tricular Electrophysiology,' Am J Physiol Heart Circ Physiol, 292:H43-H55, 2007. 

[42] T.J. Hund and Y. Rudy, 'Rate Dependence and Regulation of Action Potential and Calcium 
Transient in a Canine Cardiac Ventricular Cell Model,' Circulation 110: 3168-3174, 2004. 

[43] F. Fenton, 'Beyond Slope One: Alternans Suppression and Other Understudied Properties of 
APD Restitution,' KITP Miniprogram on Cardiac Dynamics, Kavli Institute for Theoretical 
Physics, Santa Barbara, CA, July 28, 2006. 

[44] Z. Qu and J.N. Weiss, 'The chicken or the egg? Voltage and calcium dynamics in the heart,' 
Am J Physiol Heart Circ Physiol, 293:H2054~H2055, 2007. 



19 



[45] D. Sato, Y. Shiferaw, Z. Qu, A. Garfinkel, J.N. Weiss, A. Karma, 'Inferring the cellular origin 
of voltage and calcium alternans from the spatial scales of phase reversal during discordant 
alternans,' Biophys J, 92:L33-L35, 2007. 

[46] P.N. Jordan, D.J. Christini, 'Characterizing the contribution of voltage- and calcium- 
dependent coupling to action potential stability: implications for repolarization alternans,' 
Am J Physiol Heart Circ Physiol, 293:H2109-H2118, 2007. 

[47] J.J. Fox, R.F. Gilmour Jr., and E. Bodenschatz, Conduction Block in One-Dimensional Heart 
Fibers, PRL, 89:198101-1, 2002. 

[48] J.J. Fox, M.L. Riccio, P. Drury, A. Werthman, and R.F. Gilmour Jr., 'Dynamic mechanism 
for conduction block in heart tissue,' New Journal of Physics, 5:101.1-101.14, 2003.09, 2004. 

[49] N.F. Otani, 'Theory of action potential wave block at-a-distance in the heart,' PRE, 75, 
021910, 2007. 



20 



